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Abstract: 

For a large class of quantized ergodic flows the quantum ergodicity theorem due to Shnirelman, 
Zelditch, Colin de Verdiere and others states that almost all eigenfunctions become equidistributed 
I in the semiclassical limit. In this work we first give a short introduction to the formulation of the 
quantum ergodicity theorem for general observables in terms of pseudodifferential operators and show 
that it is equivalent to the semiclassical eigenfunction hypothesis for the Wigner function in the case of 
ergodic systems. Of great importance is the rate by which the quantum mechanical expectation values 
of an observable tend to their mean value. This is studied numerically for three Euclidean billiards 
(stadium, cosine and cardioid billiard) using up to 6000 eigenfunctions. We find that in configuration 
space the rate of quantum ergodicity is strongly infiuenced by localized eigenfunctions like bouncing 
ball modes or scarred eigenfunctions. We give a detailed discussion and explanation of these effects 
using a simple but powerful model. For the rate of quantum ergodicity in momentum space we observe 
a slower decay. We also study the suitably normalized fluctuations of the expectation values around 
their mean, and find good agreement with a Gaussian distribution. 
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1 Introduction 



In quantum chaos a lot of work is devoted to the investigation of the statistics of eigenvalues 
and properties of eigenfunctions of quantum systems whose classical counterpart is chaotic. For 
ergodic systems the behavior of almost all eigenfunctions in the semiclassical limit is described 
by the quantum ergodicity theorem, which was proven in |]T], ^, ^ ^ , see also ^ for general 
introductions. Roughly speaking, it states that for almost all eigenfunctions the expectation 
values of a certain class of quantum observables tend to the mean value of the corresponding 
classical observable in the semiclassical limit. 

Another commonly used description of a quantum mechanical state is the Wigner function p|, 
which is a phase space representation of the wave function. According to the "semiclassical 
eigenfunction hypothesis" the Wigner function concentrates in the semiclassical limit on regions 
in phase space, which a generic orbit explores in the long time limit t — > oo 0, |T2|. 



For integrable systems the Wigner function W{p, q) is expected to localize on the invariant 
tori, whereas for ergodic systems the Wigner function should semiclassically condense on the 
energy surface, i.e. W{j),q) ~ vo^Ss) ^i^iP^ l) ~ -^)' where H{p,q) is the Hamilton function 
and vo\{Tie) is the volume of the energy shell defined by H{p, q) = E. 

As we will show below the quantum ergodicity theorem is equivalent to the validity of the 
semiclassical eigenfunction hypothesis for almost all eigenfunctions if the classical system is 
ergodic. Thus a weak form of the semiclassical eigenfunction hypothesis is proven for ergodic 
systems. 

For practical purposes it is important to know not only the semiclassical limit of expectation 
values or Wigner functions, but also how fast this limit is achieved, because in applications one 
usually has to deal with finite values of h, or finite energies, respectively. Thus the so-called 
rate of quantum ergodicity determines the practical applicability of the quantum ergodicity 
theorem. A number of articles have been devoted to this subject, see e.g. [|13|, 0, |15|, H, |l^ 



and references therein. The principal aim of this paper is to investigate the rate of quantum 
ergodicity numerically for different Euclidean billiards, and to compare the results with the 
existing analytical results and conjectures. A detailed numerical analysis of the rate of quantum 
ergodicity for hyperbolic surfaces and billiards can be found in [T7[ . 



Two problems arise when one wants to study the rate of quantum ergodicity numerically. First 
the fluctuations of the expectation values around their mean can be so large that it is hard or 
even impossible to infer a decay rate. This problem can be overcome by studying the cumulative 
fluctuations, 

^l(^,^) = E \{^n,A^^)-^\ , (1) 

En<E 



where {ipn, Aipn) is the expectation value of the quantum observable A, cr{A) is the mean value 
of the corresponding classical observable <j{A) and N{E) is the spectral staircase function, see 
section ^ for detailed definitions. So Si{E, A) contains all information about the rate by which 
the quantum expectation values tend to the mean value, but is a much smoother quantity than 
the sequence of differences itself. 

Secondly, since the quantum ergodicity theorem makes only a statement about almost all 
eigenfunctions (i.e. a subsequence of density one, see below), there is the possibility of not 
quantum-ergodic subsequences of eigenfunctions. Such eigenfunctions can be for example so- 
called scarred eigenfunctions ||18|, which are localized around unstable periodic orbits, or 



in billiards with two parallel walls so-called bouncing ball modes, which are localized on the 
family of bouncing ball orbits. 

Although such subsequences of exceptional eigenfunctions are of density zero, they may have a 
considerable influence on the behavior of Si{E, A). This is what we find in our numerical com- 
putations for the cosine, stadium and cardioid billiard, which are based on 2000 eigenfunctions 
for the cosine billiard and up to 6000 eigenfunctions for the stadium and cardioid billiard. 
In order to obtain a quantitative understanding of the influence of not quantum-ergodic sub- 
sequences on the rate, we develop a simple model for Si{E,A) which is tested successfully for 
the corresponding billiards. The application of this model in the case of the stadium billiard 
reveals in addition to the bouncing ball modes a subsequence of eigenfunctions, which appear 
to be not quantum-ergodic in the considered energy range. 

A further interesting question is if the boundary conditions have any influence on the rate of 
quantum ergodicity. This is indeed the case, for observables located near the boundary a strong 
influence on the behavior of Si{E,A) is observed. But for E ^ oo this influence vanishes, so 
the asymptotic rate is independent of the boundary conditions. 

After having some knowledge on the rate by which the expectation values (ipnyAipn) tend to 
their quantum-ergodic limit ^"(A) one is interested how the suitably normalized fluctuations 
{ipn.Aipn) — (^{A) are distributed. It is conjectured that they obey a Gaussian distribution, 
which we can confirm from our numerical data. 

The outline of the paper is as follows. In section |^ we first give a short introduction to the 
quantum ergodicity theorem and its implications. Then we discuss conjectures and theoretical 
arguments for the rate of quantum ergodicity given in the literature. In particular we study 
the influence of not quantum-ergodic eigenfunctions. In section ^ we give a detailed numerical 
study on the rate of quantum ergodicity for three Euclidean billiard systems for different types 
of observables, both in position and in momentum space. This includes a study of the influence 
of the boundary and a study of the fluctuations of the normalized expectation values around 
their mean. We conclude with a summary. Some of the more technical considerations using 
pseudodifferential operators are given in the appendix. 



2 Quantum ergodicity 

The classical systems under consideration are given by the free motion of a point particle inside 
a compact two-dimensional Euclidean domain i7 C with piecewise smooth boundary, where 
the particle is elastically reflected. The phase space is given by x f2, and the Hamilton 
function is (in units 2m = 1) 

Hip,q)=p\ (2) 
The trajectories of the flow generated by H{p, q) lie on surfaces of constant energy E, 

Ss:={(p,g)GR2xn \ p' = E] , (3) 

which obey the scaling property He = -E^Si := {{E^p,q) \ {p,q) E Si} since the Hamilton 
function is quadratic in p. Note that Ei is just x Q. 

The classical observables are functions on phase space x f2, and the mean value of an 
observable a{p, q) at energy E is given by 



a" = -L s / a{P, l) d/i = ^ / / Q) ^{P^ - E) dpdq , (4) 

K^e) J V0l(S£;) J J 



-E _ 

vol (Si 



where d/x = | dipdq is the Liouville measure on and vo\(T,e) = j^,^ d/x. The unusual factor 
1/2 in the Liouville measure is due to the fact that we have chosen and not as Hamilton 
function. For the mean value at energy ii^ = 1 we will for simplicity write a. 
The corresponding quantum system which we will study is given by the Schrodinger equation 
(in units h = 2m = 1) 

-Aipniq) = Eni)n{q) , q en , (5) 

with Dirichlet boundary conditions: ^niq) = for q G dil. Here A = ^ + ^ denotes the 
usual Laplacian, and we will assume that the eigenvalues are ordered a.s Ei < E2 < E3 . . . and 
that the eigenf unctions are normalized, \4'n{q)\'^ dg = 1. 

The quantum ergodicity theorem describes the behavior of expectation values {ipn, Aipn) in the 
high energy (semiclassical) limit En — > 00, and relates it to the classical mean value (^. The 
observable A is assumed to be a pseudodifferential operator, so before we state the theorem we 
have to introduce the concept of pseudodifferential operators, see e.g. PO, pTl El]. 



2.1 Weyl quantization and pseudodifferential operators 



It is well known that every continuous operator A : C^{VL) is characterized by its 

Schwarz kernel Ka € V'iVL x Vt) such that Aifji^q) = f^KA{q,q')4'{q') dq', where V'{fl) is the 
space of distributions dual to C^{fl), see e.g. chapter 5.2]. In Dirac notation one has 
KA{q,q') = {q\A\q'). With such an operator A one can associate its Weyl symbol, defined as 



W[A]ip,q):-- 



p 



2'^+ 2 ' 



(6) 



which in general is a distribution An operator A is called a pseudodifferential operator, if 
its Weyl symbol belongs to a certain class of functions. One of the simplest classes of symbols 
is ^"^(M^ X VL) which is defined as follows: a(p, q) G S'^{R^ x fi), if it is in C°^(M2 x VL) and for 
all multiindices a, j3 the estimate 

-9l"l f)\P\ m-]a| 

— ^a(;.,g) <C.A^ + \P?)^ (7) 

holds. Here m is called the order of the symbol. The main point in this definition is that 
differentiation with respect to p lowers the order of the symbol. For instance polynomials of 
degree m in ^|^,|^^ CQ,/(g)p"', whose coefficients satisfy f^CQ,/(g) < are in S''"(]R^ x 

An operator A is called a pseudodifferential operator of order m, A e S"^{Q), if its Weyl symbol 
belongs to the symbol class S""(M^ x Q), 

=^ W[v4](p,g) e 5™(M' X . 



A G s'^in) 



(8) 



For example if the Weyl symbol is a polynomial in p, then the operator is in fact a differen- 
tial operator and so pseudodifferential operators are generalizations of differential operators. 
Further examples include complex powers of the Laplacian, (-A)^/^ ^ s'^'iQ), see p5| , ^ |22[1 . 
On the other hand, to any function a G S'™'(]R" x Q) one can associate an operator Op [a] G 



Op[a]/(g) : = 



(2vr) 



P, 



q + q 



f{q') dq'dp , 



(9) 



such that its Weyl symbol is a, i.e. W[Op[a]] = a. This association of the symbol a to the 
operator A is called Weyl quantization of a. 

In practice one often encounters symbols with a special structure, namely those which have an 
asymptotic expansion in homogeneous functions in p, 



oo 



aiP,Q) ^^(^m-kiP,q), with am-kiXp,q) = X"" ''am-kiP,q) for A > . (10) 

k=0 

Note that it is not required that m is an integer, all m G M are allowed. Since the degree 



of homogeneity tends to — oo this can be seen as an expansion for \p\ oo; see pO|, ^ 
for the exact definition of this asymptotic series. These symbols are often called classical or 
polyhomogeneous, and we will consider only operators with Weyl symbols of this type. The 
space of this operators will be denoted by S^{Q). If A G S^{Q) and W{A) ~ '^'^^Qdm-k, then 
the leading term am{p,q) is called the principal symbol of A and is denoted by a{A){p,q) : = 
am{p, q)- It plays a distinguished role in the theory of pseudodifferential operators. One reason 
for this is that operations like multiplication or taking the commutator are rather complicated 
in terms of the symbol, but simple for the principal symbol. For instance one has |]20| , |2l| 

a{AB) = a{A)a{B) , a{[A, B]) = i{a{A), a{B)} , (11) 

where {■,■} is the Poisson bracket. It furthermore turns out that the principal symbol is a 
function on phase space, i.e. has the right transformation properties under coordinate transfor- 
mations, whereas the full Weyl symbol does not have this property. 

So every operator A with principal symbol cT(y4) can be seen as a quantization of the classical 
observable ct(A). The existence of different operators with the same principal symbol just 
reflects the fact that the quantization process is not unique. Furthermore, one can show that 
the leading asymptotic behavior of expectation values of such operators for high energies only 
depends on the principal symbol, as it should be according to the correspondence principle. 
This is a special case of the Szego limit theorem, see |^ chapter 29.1]. 

One advantage of the Weyl quantization over other quantization procedures is that the Wigner 
function of a state \ip) appears naturally as the Weyl symbol of the corresponding projection 
operator 

W [\^) m {p, q) = j e^'^'P ^(^q-l^^p(^q+^-^ dq' . (12) 

In the following we will use for a Wigner function of an eigenstate ipn the simpler notation 
Wn(p, g) := W [!'?/'„)('?/'„ I] (p, g). For the expectation value {ip,Aip) one has the well known 
expression in terms of the Weyl symbol W[A] and the Wigner function W [|^)('?/'|], 

(^,A^) = ^ jj W[A]ip,q)W[\^){i;\]ip,q)dpdq . (13) 

Pseudodifferential operators of order zero have a bounded Wigner function, and therefore a 
bounded principal symbol cr{A); this boundedness of the classical observable carries over to the 
operator level: the operators in S^{Q) are bounded in the L^-norm. 

The definition of pseudodifferential operators can be generalized to manifolds of arbitrary di- 
mension, the preceeding formulas are then valid in local coordinates. The symbols of these 



operators only live in local charts, but the principal symbols can be glued together to a func- 
tion on the cotangent bundle T*Q which is the classical phase spaceQ. 



2.2 Quantum limits and the quantum ergodicity theorem 

In quantum mechanics the states are elements of a Hilbert space, or more generally linear 
functionals on the algebra of observables. In classical mechanics the pure states are points 
in phase space, and the observables are functions on phase space. More generally the states 
are measures on phase space, which are linear functionals on the algebra of observables. The 
pure states are then represented as delta functions. The eigenstates of a Hamilton operator are 
those which are invariant under the time evolution defined by H. In the semiclassical limit they 
should somehow converge to measures on phase space which are invariant under the classical 
Hamiltonian flow. The measures which can be obtained as semiclassical limits of quantum 
eigenstates are called quantum limits. 

More concretely the quantum limits can be described as limits of sequences of Wigner functions. 
Let {V'njneN be an orthonormal basis of eigenf unctions of the Dirichlet Laplacian —A, and 
{W„}„gN the corresponding set of Wigner functions, see equation ([l2|) . We first consider 
expectation values for operators of order zero, and then extend the results to operators of 
arbitrary order. 

Because pseudodifferential operators of order zero are bounded, the sequence of expectation 
values {{ipn, Aipn)}nm is bounded too. Every function a G C°°(Si) can be extended to a 
function in C°°(M^\{0} x by requiring it to be homogeneous of degree zero in p. Via the 
quantization Op[a]0 of a and equation ( |I3| ) one can view the Wigner function W„(p, g) as a 
distribution on C°°(Ei), 

a {ijjn, Op[a]^pn) = jj a{p, q) Wn{p, q) dpdg . (14) 

The sequence of these distributions is bounded because the operators Op [a] are bounded. The 
accumulation points of {W„(p, g)}neN are called quantum limits //^(p, Q'), and we label them by 
k E I, where / is some index-set. Corresponding to the accumulation points fik{p,(l), the se- 
quence {Wn{p,q)}nm can be split into disjoint convergent subsequences IJfce/'{^ri'=(P5 Q')}ieN = 
{W„(j9, g)}„gN- That is, for every k we have 

jj a{p, q) W„fc (p, q)dpdq = jj a{p, g)/ifc(p, q) dp dq , (15) 

nxIR2 OxIR2 

for all a G C°°(Ei) viewed as homogeneous functions of degree zero on phase space. This 
splitting is unique up to a finite number of terms, in the sense that for two different splittings 



^If one wants to realize the semiclassical limit not as the high energy limit, but as the limit of ?i — s- one 
has to incorporate H explicitly in the quantization procedure. In the framework of pseudodifferential operators 
this has been done by Voros in |l^ , see also ^ . 

^Strictly speaking is a not an allowed Symbol because it is not smooth a.t p — 0. Let x{p) G C°°(M^) satisfy 
xip) = for \p\ < 1/4 and xip) — 1 for \p\ ^ 1/2- By multiplying a with this excision function xip) we get 
a symbol xi S ^^(M^ x il), whose Weyl quantization Op[xa] is in ^^(rj). But the semiclassical properties of 
Op[xa] are independent of the special choice of xip)i which can be see e.g. in (p^, since Wn is concentrated on 
the energy shell T,e„ for n — > oo. Therefore we will proceed for simplicity with a instead of xa. 



the subsequences belonging to the same accumulation point differ only by a finite number of 
terms. 

As has been shown in the quantum limits fik are measures on which are invariant 



under the classical flow generated by H{p, q). 

One of the main questions in the field of quantum chaos is which classical invariant measures 
on C°°(Si) can actually occur as quantum limits of Wigner functions. E.g., if the orbital 
measure along an unstable periodic orbit occurs as quantum limit fi^, then the corresponding 
subsequence of eigenfunctions has to show an enhanced probability, i.e. scarring, along that 
orbit. 

Given any quantum limit fik one is furthermore interested in the counting function Nk{E) : = 
^{E^k < E} for the corresponding subsequence {W„fc}jgN of Wigner functions. Since the 
subsequence {W„fc}jgN is unique up to a finite number of elements, the corresponding counting 
function Nk{E) is unique up to a constant. 

One should keep in mind that we have defined the quantum limits and their counting functions 
here with respect to one chosen orthonormal basis of eigenfunctions {ipn{Q)}neN- If one takes a 
different orthonormal base of eigenfunctions {^n(Q')}neN, the counting functions corresponding 
to the quantum limits, or even the quantum limits themselves, may change. So when studying 
the set of all quantum limits one has to take all bases of eigenfunctions into account. 
The lift of any quantum limit from Si to the whole phase space x follows straightforward 
from some well known methods in pseudodifferential operator theory, as shown in appendix 0. 
For a pseudodifferential operator of order m, A E S'^{Q), one gets for the expectation values 

lim E~^^'^{ip^k,Aip^k) = /ifc((T(A)|si) = / (T(A)(p,g)/ifc(p,g)d/i . (16) 

Si 

In terms of the Wigner functions this expression can be written as (see appendix |C[) 

lim 4 W„.(4p, Q) = f^kip, q) ~ . (17) 

Without the scaling of p with \/E we have 

5{H{p,q)-E^k) 

W„fe q) ~ /ifc p, q) — , (18 

vol(Si5^J 

for E^k — >■ oo, and fik{p, <?) is extended from Si to the whole phase space by requiring it to be 
homogeneous of degree zero in p. 

For ergodic systems the only invariant measure whose support has nonzero Liouville measure 
is the Liouville measure itself. For these systems the quantum ergodicity theorem states that 
almost all eigenfunctions have the Liouville density as quantum limit. 
Quantum ergodicity theorem 



Let C be a compact 2- dimensional domain with piecewise smooth boundary, and let {tpn} 
be an orthonormal set of eigenfunctions of the Dirichlet Laplacian A on Q. If the classical 
billiard flow on the energy shell Si = S*^ x is ergodic, then there is a subsequence {nj} C N 
of density one such that 



lim (V'n,.,AV'„,) =aiA) , (19) 



for every polyhomogeneous pseudodijferential operator A G S^{fl) of order zero, whose Schwarz 
kernel Kyi{q,q') = {q\A\q') has support in the interior of Q x Q. Here cr{A) is the principal 
symbol of A and <j{A) is its classical expectation value, see eq. 
A subsequence {uj} C N has density one if 

lim ' ' "\ = 1 , (20 

s-oo N{E) ' ^ ^ 

where N{E) := | En < E} is the spectral staircase function, counting the number of 

energy levels below a given energy E. So almost all expectation values of a quantum observable 
tend to the mean value of the corresponding classical observable. 

The special situation that there is only one quantum limit, i.e. the Liouville measure, is called 
unique quantum ergodicity. This behavior is conjectured to be true for the eigenfunctions of 
the Laplacian on a compact manifold of negative curvature |T5|. 

We have stated here for simplicity the quantum ergodicity theorem only for two dimensional 
Euclidean domains, but it is true in far more general situations. For compact Riemannian 
manifolds without boundary the quantum ergodicity theorem was given by Shnirelman 
Zelditch and Colin de Verdiere [Q. For a certain class of manifolds with boundary it was 
proven in [^, without the restriction on the support of the Schwarz kernel of the operator 
A. The techniques of can possibly be used to remove these restrictions here as well, see 
the remarks in |^0[. One can allow as well more general Hamilton operators; on manifolds 
without boundary every elliptic selfadjoint operator in S'^i{Q) is allowed, and on manifolds 
with boundary at least every second order elliptic selfadjoint differential operator with smooth 
coefficients is allowed. This includes for instance a free particle in a smooth potential or in a 
magnetic field. In the semiclassical setting, where the Hamilton operator and the observables 
depend explicitly on h, a similar theorem for the limit h —>■ has been proven in ^j, see also 
1^ for an introduction. 

In light of the correspondence principle the quantum ergodicity theorem appears very natural: 
Classical ergodicity means that for a particle moving along a generic trajectory with energy 
E, the probability of finding it in a certain region [/ C S^; of phase space is proportional to 
the volume vol(?7) of that region, but does not depend on the shape or location of U. The 
corresponding quantum observable is the quantization of the characteristic function xu of U, 
and by the correspondence principle one expects that the expectation value of this observable 
in the state ipn tends to the classical expectation value for En oo. And this is the content of 
the quantum ergodicity theorem. 

In terms of the Wigner functions W„ the theorem gives, see eq. ([T^), 

6iHip,q)-E^^) 

(p, q) — , (21 

vo1(Sb„^) 

for j — >• oo, for a subsequence {uj} C N of density one. So almost all Wigner functions 
become equidistributed on the energy shells T,e„.- That is, for ergodic systems the validity of 
the semiclassical eigenfunction hypothesis for a subsequence of density one is equivalent to the 
quantum ergodicity theorem. 



2.3 Some examples 



As an illustration of the quantum ergodicity theorem, and for later use, we now consider some 
special observables whose symbol only depends on the position q or on the momentum p. 



If the symbol only depends on the position q, i.e. a{p,q) = a{q), the operator is just the 
multiplication operator with the function a{q), and one has 



(V',A^/.) = {tp,aij) = / a{q) \ij{q)\' dq . (22) 



In the special case that one wants to measure the probability of the particle to be in a given 
domain D C Q, the symbol is the characteristic function of D, i.e. a{p,q) = Xd{q)- Then 
Op[Xd] is not a pseudodifferential operator, but nevertheless the quantum ergodicity theorem 
remains valid for this observable [Q. Since the principal symbol is then cr{A) = xd we obtain 
for its mean value 

Thus the quantum ergodicity theorem gives for this case 

hm / \^^^(q)\^dq = ^^^^ (24) 



D 



vol (n) 



for a subsequence {rij} C N of density one. As discussed at the end of the previous section this 
is what one should expect from the correspondence principle. 

If instead the symbol depends only on the momentum p, i.e. a{p,q) = a{p), one obtains from 
for the expectation value 

(^,AV')= / a{p)mp)\'dp . (25) 



In the same way as in ^ for a characteristic function in position space, it follows that the quan- 
tum ergodicity theorem remains valid for the case where a{p) = Xc{eA9)ip) is the characteristic 
function of a circular sector in momentum space of angle 6. In polar coordinates this is given 
by the set 

C{e,Ae):= {{r,^) \ r eR+, ^e[e-Ae/2,e + Ae/2]} . (26) 
The mean value of the principal symbol then reduces to 

"^^^ = J ^^(^.^^)(^) = ^ ' (27) 

which does not depend on 9. Thus the quantum ergodicity theorem reads in the case of a 
characteristic function in momentum space 

//\f) 
fc,(p)pdp=— (28) 

c(e,Af) 



for a subsequence {uj} C N of density one. This means that quantum ergodicity implies an 
asymptotic equidistribution of the momentum directions of the particle. 



It is instructive to compute the observables discussed above for certain integrable systems. 
First consider a two-dimensional torus. The eigenfunctions, labeled by the two quantum num- 
bers n, m G Z, read ipn,m{x,y) = exp(27rinx) exp(27rim?/). Obviously, these are "quantum- 
ergodic" in position space, since \ilJn,m{x,y)\'^ = 1, but they are not quantum-ergodic in mo- 
mentum space. Even in position space the situation changes if one takes a different orthogonal 
basis of eigenfunctions (note that the multiplicities tend to infinity), see for a complete 
discussion of the quantum limits on tori. A similar example is provided by the Dirichlet or 
Neumann eigenfunctions of a rectangular billiard. 

The circle billiard shows a converse behavior. Let the radius be one, then the eigenfunctions 
are given in polar coordinates by 

^M{r,(l)) = NkiJi{jk,ir)e''^ . (29) 

Here jk,i is the fc-th zero of the Bessel function Ji{x), x > 0, and Nki is a normalization constant. 
These eigenfunctions do not exhibit quantum ergodicity in position space. But for their Fourier 
transforms one can show that 



dp = ^ , (30) 



C(6»,A6») 



and so we have "quantum ergodicity" in momentum space. 

A remarkable example was discussed by Zelditch |Q. He considered the Laplacian on the 
sphere S"^. Since the multiplicity of the eigenvalue /(/ + 1) is 2/ + 1, which tends to infinity as 
Z — > oo, one can choose infinitely many orthonormal bases of eigenfunctions. Zelditch showed 
that almost all of these bases exhibit quantum ergodicity in the whole phase space. Although 
this is clearly an exceptional case due to the high multiplicities, it shows that one has to be 



careful with the notion of quantum ergodicity. In a recent work Jakobson and Zelditch 
have furthermore shown that for the sphere all invariant measures on phase space do occur 
as quantum limits. One might conjecture that for an integrable system all classical measures 
which are invariant under the flow and all symmetries of the Hamilton function do occur as 
quantum limits. 

The general question whether quantum ergodicity for all orthonormal bases of eigenfunctions 
in the whole phase space implies ergodicity of the classical system is still open. 



2.4 The rate of quantum ergodicity 

We now come to the central question of the approach to the quantum-ergodic limit. First we 
note that an equivalent formulation of the quantum ergodicity theorem, which avoids choosing 
subsequences, is given by 



lim — ; — - 

E^oc N{E) 



J2 \{^n, A^n) - (r{A) 

E„<E 







(31) 



This equivalence follows from a standard lemma concerning the influence of subsequences of 

density zero on the average of a sequence, see e.g. [^, Theorem 1.20]. 

In order to characterize the rate of approach to the ergodic limit the quantities 



SUE, A) 



N{E) 



Utfjn^A^.,) -a{A) 



E„<E 



(32) 



have been proposed and studied in |]T3|, |T^. Quantum ergodicity is equivalent to Sm{E,A) — > 
for — >■ oo and m > 1. 

Let us first summarize some of the known resuhs for the rate of quantum ergodicity. 
Zelditch proved in |]l3l by relating the rate of quantum ergodicity to the rate of convergence 
of classical expectation values, and using a central limit theorem for the classical flow, that for 
compact manifolds of negative curvature SmiE,A) = 0((logi?)~™/^). However this bound is 
believed to be far from being sharp. Moreover in lower bounds for Sm{E,A) have been 
derived. In |3^, |3^ it is proven for a Hecke basis of eigenfunctions on the modular surface 
that S2{E,A) < C{e)E~2^'^ for every e > 0. It is furthermore conjectured 0, |I^ that this 



estimate is also valid for the eigenfunctions of the Laplacian on a compact manifolds of negative 
curvature, and moreover that it is satisfied for each eigenstate individually: | {ipn, Aipn) — cr(^) I < 
C(e)-E~4+^ for every e > 0. 

In a study of 5*2 A) based on the Gutzwiller trace formula has been performed. For 
completely desymmetrized systems having only isolated and unstable periodic orbits, the so- 
called diagonal approximation for a double sum over periodic orbits and further assumptions 
lead to 

Here g = 2 if the system is invariant under time reversal, and otherwise g = 1, and p{A) is the 
variance of the fluctuations of = ^ J^'' a{A){'y{t)) dt around their mean <y{A), computed 

using all periodic orbits 7 of the system. More precisely, it is assumed that \A^ — cr(A)P ~ 
p{A)/T^, where denotes the primitive length of 7. 

In the general case where not all periodic orbits are isolated and unstable it is argued that the 
rate of quantum ergodicity is related to the decay rate of the classical autocorrelation function 
C(r) [g. If C(r) ~ T-'i then the result is 

Th I for 77 > 1 

S2iE, A)^ j C(r) dr ~ J In i^-^E^I^'^E-^'^ for = 1 , (34) 
forr/<l 

where = jg ^j-^g so-called Heisenberg time. 

For the stadium billiard and the Sinai billiard [39| it is believed that the correlations decay 



as ~ 1/r, see and |^ for numerical results for the Sinai billiard. Thus for both the stadium 
and the Sinai billiard a logarithmic contribution to the decay of 6*2 (-E, A) is expected. 
Also a Gaussian random behavior of the eigenfunctions [|Tl|] implies in position space a rate 
S2{E,A) = 0{E-2), which follows from [||, chapter IV], see also [^, 



Random matrix theory (see section VII]) predicts for suitable observables the same rate 
S2{E, A) = 0{E~2^^ and furthermore Gaussian fluctuations of {{ipn, Aijjn) ~'^{^))/ a/ S2{En, A) 
around zero, which we study numerically in section |3.3| . 

Since for the systems under investigation we have not quantum-ergodic subsequences of eigen- 
functions, we now discuss in general the influence of such subsequences on the behavior of 
Si{E, A). To this end we split the sequence of eigenfunctions into two subsequences. The first, 
denoted by {i'n'}, contains all quantum-ergodic eigenfunctions, i.e. the corresponding quantum 
limit of the associated sequence of Wigner functions is the Liouville measure. The counting 
function of this subsequence will be denoted by N'{E). The other sequence {ipn"} contains all 



not quantum-ergodic eigenfunctions. This subsequence may have different quantum limits fik 
which are all different from the Liouville measure. Their counting function will be denoted by 
N"{E). Examples would be a subsequence of bouncing ball modes or eigenfunctions scarred 
by an unstable periodic orbit. Similarly we split Si{E,A) into two parts corresponding to the 
two classes of eigenfunctions. Due to the separation N{E) = N'{E) + N"{E) we obtain 



Si{E,A) 



N{E) 



En<E 



^'(^\s',iE,A) + ^S';iE,A) 



NiE) 



NiE) 



l-^)s[iE,A) + ^S';iE,A). 



N{E) 



NiE) 



Here we defined 



S[{E,A) 
S'liE,A) 



N'{E) 



E„,<E 



N"{E) 



E,,<E 



(35) 

(36) 
(37) 



So the behavior of Si{E,A) is given in terms of the three quantities S[{E,A), S'l^E^A) and 
N"{E), which describe the behavior of the quantum-ergodic and the not quantum-ergodic 
subsequences, respectively. 

The behavior of Si{E, A) can be described in terms of the non ergodic quantum limits and 
their counting functions. We split the not quantum-ergodic subsequence into convergent sub- 
sequences corresponding to the quantum limits fik ^ /i, {ipn"} = Ufc{V'n'=}i6N, with N"{E) = 



Nk{E), and {ipn'^, Ai/j^k) —a{A) ~ /x^ ya{A) — a{A) ) . Then S'{{E, A) is asymptotically given 

by 

1 



S'I{E,A) 



J2NkiE)\f-ik{cxiA)-a{A) 



and the limit 



u"{A) := lim S';{E,A) 

E^oo 



(38) 



(39) 



only depends on cr{A) and defines an invariant measure on Si. 

Let us assume for the quantum-ergodic part of Si{E, A) a certain rate of decay, 

S[iE, A) = u'{A)E-'' + o{E~'') , 

and for the counting function of the not quantum-ergodic states 

N"{E) = cE^ + o{E^) , 



(40) 



where by quantum ergodicity a > and (3 < 1. With Weyl's law N{E) 



then obtain in eq. 



for Si{E,A) 



vol(Q) 
4tt 



(41) 

E + 0{E^) we 



Si{E,A) = u'{A)E~'' + 



Attc 
vol(fi) 



"{A)E^-^ + o{E~'') + o{E^"^) . 



(42) 



One sees that if —a > /3 — 1, the asymptotic behavior of S\{E^ A) is governed by the quantum- 
ergodic sequences of eigenf unctions, whereas in the opposite case, —a </? — !, the not quantum- 
ergodic sequences dominate the behavior asymptotically. Especially if/3 — 1 > —1/4, i.e. 
fi > 3/4, the rate of quantum ergodicity cannot be 0{E~^). 

To obtain a simple model for the rate of quantum ergodicity, let us now assume that the 
conjectured optimal rate is valid for the subsequence of quantum-ergodic eigenfunctions, that 
is a = 1/4 can be chosen in eq. (^Ol). To be more precise it should be S[{E, A) = 0{E~^^^^^) for 
every e > 0, but for comparison with numerical data on a finite energy range we will assume that 
e = 0. For the not quantum-ergodic eigenfunctions the knowledge of their counting function 
N"{E) is very poor; in general it is unknown. Thus if we neglect the higher order terms in 
eqs. (|40|) and ( Pl) we obtain from ( ^5l) and (|39|) a simple model for the behavior of Si{E, A), 

sr"(E, A) = (i - A-m-i + ^-'V)^-' . (43) 

The first factor in braces will only be important if jS is close to 1. 

We will now discuss the influence of a special type of not quantum-ergodic subsequences in 
more detail. In billiards with two-parallel walls, one has a subsequence of so-called bouncing 
ball modes [45|, which are localized on the bouncing ball orbits, see fig. |l]b) for an example of 



such an eigenf unction. In our previous work we showed that for every /5 < 1 there exists 
an ergodic billiard which possesses a not quantum-ergodic subsequence, given by bouncing ball 
modes, whose counting function is asymptotically of order E^. But for (3 = 1 — 6, with some 
small S > 0, equation (|^) shows that Si{E,A) = 0{E~^) at least for some A. So the best 
possible estimate on the rate of quantum ergodicity which is valid without further assumptions 
on the system other than ergodicity is 

Si{E,A) = o(l) , i.e. lim Si{E,A) = . (44) 

E^oo 

Especially for the Sinai billiard the result for the exponent is /5 = 9/10 and therefore Si{E, A) ~ 
which contradicts the result (|33) from ||16||. 



If the bouncing ball modes are the only not quantum-ergodic eigenfunctions, or at least consti- 
tute the dominant contribution to them, then N'\E) ~ A^bb(-E') ~ cE^ . The exponent /3 and 



jy"{A) are explicitly known, and the constant c is known from a numerical fit in for the 



billiards we will consider in the next section. Thus in this case the only free parameter in the 



model (H) is u'{A). 

The asymptotic behavior of ( ^Hj ) is governed by the term with the larger exponent, but this 
can be hidden at low energies if one of the constants is much larger than the other. Assume 
for instance that {3 — 1 > 1/4, i.e. the not quantum-ergodic eigenfunctions dominate the rate 
asymptotically. If 

4W(A) ^ ,^ , 

^ ^ < 1 , (45) 



vo\{ny{A) 



for an observable A, then up to a certain energy Si{E,A) will be approximately proportional 
to E~^. In numerical studies where only a finite energy range is accessible such a behavior 
can hide the true rate of quantum ergodicity. This will be seen most drastically for the cosine 
billiard, see section |3.1.1| . 

The main ingredient of the model (|i3| ) is the conjectured behavior of the rate for the quantum- 
ergodic eigenfunctions. By comparing (^31) with numerical data for different observables one 



can test this conjecture. If this conjecture is true then it means that the only deviations 
from the optimal rate of quantum ergodicity are due to subsequences of not quantum-ergodic 
eigenfunctions. 

Clearly similar models based on a splitting like (|35|) can be developed for other situations as 
well. E.g., if the eigenfunctions split into a quantum-ergodic subsequence of density one with 
rate proportional to E""^/^ and a quantum-ergodic subsequence of density zero with a slower, 
and maybe spatial inhomogeneous, rate, one would expect a similar behavior of 5*1 [E, A) as 
in the case considered above. So it will be hard without some a priori information on not 
quantum-ergodic eigenfunctions to distinguish between these two scenarios. 



3 Numerical results 

In order to study the rate of quantum ergodicity numerically we have chosen three different 
Euclidean billiard systems, given by the free motion of a point particle inside a compact domain 
with elastic reflections at the boundaries. See fig. |^ for the chosen billiard shapes. 



The first is the stadium billiard, which is proven to be ergodic, mixing and a i^-system [pq . 
The height of the desymmetrized billiard is chosen to be 1, and a denotes the length 
of the upper horizontal line. For this system our analysis is based on computations of the 
first 6000 eigenfunctions for odd-odd parity, i.e. everywhere Dirichlet boundary conditions in 
the desymmetrized system with parameter a = 1.8. We also studied stadium billiards with 
parameters a = 0.5 and a = 4.0 using the first 2000 eigenfunctions in each case to investigate 
the dependence on a, see below. The stadium billiard is one of the most intensively studied 



systems in quantum chaos, for investigations of the eigenfunctions see e.g. |^ |18|, |T9|, |S 
and references therein. 

The second system is the cosine billiard, which is constructed by replacing one side of a rect- 
angular box by a cosine curve. The cosine billiard has been introduced and studied in detail in 



50| , |51[| . The ergodic properties are unknown, but numerical studies do not reveal any stability 
islands. If there were any they are so small that one expects that they do not have any influence 
in the energy range under consideration. The height of the cosine billiard is 1 and the upper 
horizontal line has length 2 in our numerical computations. The cosine is parameterized by 
B{y) = 2 + ^(1 + cos(7r?/)), see fig. ^). For our analysis of this system we used the first 2000 
eigenfunctions with Dirichlet boundary conditions everywhere. 



Figure 1: Left: Density plots \ipn{(l)\'^ for three different odd-odd eigenfunctions 
of the a = 1.8 stadium billiard: a) n = 1992, "generic" h) n = 1660, bouncing 
ball mode c) n = 1771 localized eigenfunction. Right: Density plots for two eigen- 
functions of the cardioid billiard with odd symmetry: d) n = 1816, "generic" e) 
n = 1817, localized along the AB orbit. Notice that according to the quantum 
ergodicity theorem the non-localized eigenfunctions of type a) and d) are the over- 
whelming majority. 



The third system is the cardioid billiard, which is the limiting case of a family of billiards 
introduced in [^. The cardioid billiard is proven to be ergodic, mixing, a ii'-system and a 

Both the classical system ||52|, ^ Q and the quantum 
I, ^ . The eigenvalues of the cardioid 



Bernoulli system [53, 54, 55, 56, 57 



mechanical system have been studied in detail 1^1], ^ 



billiard have been provided by Prosen and Robnik and were calculated by means of the 



conformal mapping technique, see e.g. pT], BSl BH]. Using these eigenvalues, our study is based 



on computations for the first 6000 eigenfunctions of odd symmetry, which were obtained from 
the eigenvalues by means of the boundary integral method ^ using the singular value 
decomposition method [^. The boundary integral method was also used for the computations 
of the eigenvalues and eigenfunctions of the stadium and the cosine billiard. 
Let us first illustrate the structure of wave functions by showing density plots of \il>n{<l)\'^ for 
three different types of wave functions of the stadium billiard and two different types of the 
cardioid billiard. Fig. ^) shows a "generic" wave function, whose density looks irregular. 
Example b) belongs to the class of bouncing ball modes, and its Wigner function is localized 
in phase space on the bouncing ball orbits, see the discussion in section |2.3| . Fig. |l|c) is another 
example of an eigenfunction showing some kind of localization. Fig. |l|d) shows a "generic" 
wave function for the cardioid billiard and [IJe) is an example of an eigenfunction, which shows a 
strong localization in the surrounding of the shortest periodic orbit (with code AB, see 



We should emphasize that according to the quantum ergodicity theorem the overwhelming 
majority of states in the semiclassical limit are of the type a) and d), which we also observe for 
the eigenfunctions of the studied systems. 



3.1 Quantum ergodicity in coordinate space 

The quantum ergodicity theorem applied to the observable with symbol a{q) 



in section 2.3, states that the difference 



di{n) 



vol (A) 
vol (O) 



Xoil), discussed 



(46) 



vanishes for a subsequence of density one. The first set of domains Di for which we investigate 
the approach to the ergodic limit is shown in fig. |^. Plots of di{n) for domain D4 of the stadium 
billiard and D5 of the cardioid billiard in fig. ^ show quite large fluctuations around zero. In 
particular for the stadium billiard there are many states for which di{n) is quite large and ^4(77.) 
is quite small. As one would expect, a large number of them are bouncing ball modes. The 
fluctuations of di{n) for the cosine billiard behave similarly to the stadium billiard. 
When trying to study the rate of the approach to the quantum-ergodic limit numerically one 
therefore is faced with two problems. On the one hand di{n) is strongly fluctuating, which 
makes an estimate of the approach to the mean very difficult, if not impossible for the available 
numerical data. On the other hand one does not know a priori which subsequences should be 
excluded in (|46|). Therefore the investigation of the asymptotic behavior of the "cumulative" 
version (|3TD of the quantum ergodicity theorem is much more appropriate. For the observable 
Xoiq) we have 



1 



NiE) 



E„<E 



vol (D) 



vol (fi) 



(47) 



Figure 2: Shapes of the billiards studied numerically in this work: a) desymmetrized 
stadium billiard, b) desymmetrized cosine billiard and c) desymmetrized cardioid 
billiard. The rectangles in the interior of the billiards mark the domains Di of 
integration for studying the rate of quantum ergodicity in configuration space. 



Figure 3: Plot of di{n) = p dg— for domain 4 in the stadium billiard 

and for domain 5 in the cardioid billiard. Since |'0n(g')P > one has di{n) > — ■ 
For domain in the stadium this lower bound is attained by the bouncing ball 
modes whose probability density |'i/'n(9)P nearly vanishes in D4; they are responsible 
for the sharp edge seen in the plot of d^ln). 

In figs. ^, 1^ and [TT] we display Si{E, Xd^) for the different domains Dj, shown in fig. ^, in the 
desymmetrized cosine, stadium and cardioid billiard, respectively. One nicely sees that the 
numerically determined curves for Si{E, XOi) decrease with increasing energy. This is of course 
expected from the quantum-ergodic theorem, however since this is an asymptotic statement, it 
is not clear a priori, whether one can observe such a behavior also at low energies. It should be 
emphasized that fig. ^ is based on the expectation values {ipn, XOii^n) for 2000 eigenfunctions 
and figs. ^ and |lT] are based on 6000 eigenfunctions in each case. 
In order to study the rate of quantum ergodicity quantitatively a fit of the function 

Sl\E) = aE-i/^+^ (48) 

to the numerical data for Si{E,XDi) is performed. As discussed in section B]^, for certain 



systems a behavior Si{E, A) = 0{E~^^^'^^) for all e > is expected, so that the fit parameter e 
characterizes the rate of quantum ergodicity. A positive value of e thus means a slower decrease 
of Si{E, A) than the expected E~^^'^. The results for e are shown in tables and the insets 
in figs, ^-pl show the same curves Si{E, xOi) in a double-logarithmic plot together with these 
fit curves. The agreement of the fits with the computed functions Si{E,XDi) is very good. 
However, e is not small for all domains Di of the considered systems, rather we find several 
significant exceptions, which will be explained in the following discussion. 

3.1.1 Cosine billiard 

For the cosine billiard one would expect a strong influence of the bouncing ball modes on the 
rate, since their number increases according to as iVbb(-E) ~ c E^^^^. But the prefactor c 



turns out to be very small and therefore the influence of the bouncing ball modes is suppressed 
at low energies. The model for Si{E,A), equation (|^), gives for the cosine billiard 

Sr^'^'iE, xdJ = (1 - 0.201 u'{xD,)E-i + 0.201 uUxdJE~'-'' , (49) 



where we have inserted the values c = 0.04 and (3 = 0.87, obtained in |^ from a fit to 
A^bb(-£') which was performed over the same energy range which we consider here. For sake 
of completeness we have included the first factor, (1 — 0.201 E~^'^^), but the numerical fits we 
perform below only change marginally if one sets this factor equal to 1. 



Figure 4: Plot of Si{E,XDi) for different domains Di for the cosine billiard using 
the first 2000 eigenfunctions, see fig. ^) for the location of the domains Di. The 
inset shows the same curves in double-logarithmic representation together with a 
fit of Sf\E) = aE~^/^+^ to the numerical data. 



domain 


rel. area 


e 


a 






1 


0.018 


-0.002 


0.052 


0.0525 


0.0045 


2 


0.018 


+0.012 


0.026 


0.0468 


0.0067 


3 


0.008 


+0.013 


0.043 


0.0297 


0.0020 


4 


0.008 


+0.022 


0.023 


0.0273 


0.0030 


5 


0.015 


+0.020 


0.050 


0.0543 


0.0150 


6 


0.336 


+0.009 


0.258 


0.2471 


0.0840 


7 


0.512 


+0.023 


0.352 


0.2920 


0.1280 


8 


0.648 


+0.009 


0.381 


0.3410 


0.1620 


9 


0.800 


+0.054 


0.279 


0.3264 


0.2500 



Table 1: Rate of quantum ergodicity for the cosine billiard with domains Di as 
shown in fig. ^ and fig. |^ and in the inset of fig. |^. Shown are the results for e and 
a of the fit of Sf^{E) = aE~^^^^^ to the numerical data. Also tabulated are the 
values for the relative area of the corresponding domains, the quantities t'bb(XA) 
computed according to (^Tj) and the result i^'{xDi) of the fit of the model ( ^91) to 



The asymptotic behavior of the probability density \'4'n"{(l)\'^ of the bouncing ball modes is (in 
the weak sense) 

2 Jl/vol(-R) for g G -R „ , , 

10 for g e \l\R 

where R denotes the rectangular part of the billiard. So the expectation values are asymptoti- 
cally {ipn",XDipn") ~ vo\{D fl R) / vol(i?), and since t'bb(XD) = ^^^e^oo S" {E , xd) is the mean 
value of \{ipn", Xni^n") — yo\{D)/ vol(fi)| over all bouncing ball modes one has 



KhiXo) 



vol{D n R) vo\{D) 



(51) 



vol(i?) vol(fi) 

For fixed volume vol{D) the quantity t'bb(X-D) is maximal for domains D lying entirely outside 

vol(D) 



of the rectangular region, fbb(Xi?) = ^°\\n) ■ For domains lying entirely inside the rectangular 



part of the billiard, we have the minimal value t'bb(X-D) = i Toi(n] ■ Therefore the strongest 
contribution of the bouncing ball modes to Si{E,xd) in eq. (|^) is expected for the domains 
outside the rectangular region. 

The values for t'bb(X-Dj given in table H The largest values for the small domains are ob- 
tained for the domains outside the rectangular part of the billiard for which also the rate of quan- 
tum ergodicity is the slowest. Furthermore we see from table H that the factor 0.201 t'bb(X-Di) in 
front of E~^'^^ in equation ( ^9|) is for all domains much smaller than the prefactor a from the fit 
to (|48|). This already indicates that the contribution of the bouncing ball modes is suppressed, 
explaining why the rate for the cosine billiard is in such a good agreement with e = 0. 



Figure 5: Plot of Si{E,XDi) for two further domains Ds and Dg (dashed curve) 
in the cosine biUiard using the first 2000 eigenfunctions. Also shown is the fit 



In order to test this quantitatively we have performed a fit of the model (^) to the numerical 
data, where the only free parameter is z/'(xdJ- The accuracy of the fits is very good and the 
results for v'^XOi) are shown in table they are much larger than the corresponding prefactors 
0.201 t'bb(X-Di) of the bouncing ball part of Si{E, XdJ- Therefore the influence of the bouncing 
ball modes on the rate is negligible small on the present energy interval, despite the fact that 
asymptotically they should dominate the rate. 

The domains C Di and C D2 show a slightly slower rate than Di and D2, respec- 
tively. This is due to the fact that choosing a smaller domain D implies larger fluctuations of 
{4'n,XD4'n) for the same set of eigenfunctions. 

As an additional test we have computed Si{E, xdJ numerically for four further domains (shown 
in the inset of fig. ^), having a much larger area than the previous ones. For these domains 
^hhiXDi) is larger, and one therefore expects a stronger influence of the bouncing ball modes 
and correspondingly a slower rate of quantum ergodicity. The results are shown in table ^ 
and fig. ^ and our findings are completely consistent with the previous ones as well as with 
the model (^). We also observe in fig. ^ that for the large domains, except for the whole 
rectangular part Dg = R, the rate is faster at low energies, than at high energies. This is due 



to the influence of the boundary and will be discussed in section p. 1.4 . 

Summarizing the results for the cosine billiard, we found that the rate of quantum ergodicity 
is in excellent agreement with a rate proportional to for the subsequence of quantum- 

ergodic eigenfunctions. The phenomenological model S^°'^'^^{E,xd), eq. (^9]), is in very good 
agreement with the numerical data, especially in view of the fact that it contains only one free 
parameter. Furthermore, the cosine billiard provides an impressive example of a system for 
which the asymptotic regime for 5*1 [E, A) is reached very late. Up to the 2000th eigenfunction 
the asymptotic behavior Si{E,A) ~ CE~^^^^ is almost completely hidden. A continuation of 
^modei^^^ domain R = Dg with the strongest influence of the bouncing ball modes, 

shows that at ~ 10® the two contributions have the same magnitude, and one has to go up 
as high as ~ 10^° to see the asymptotic behavior Si{E,xr) ~ CE~^^^^. Therefore there is 
no contradiction between the observed fast rate of quantum ergodicity in the present energy 
range and the increase of the number of bouncing ball modes iVbb(-E) ~ c E^^^^ found in . 



3.1.2 Stadium billiard 

For the stadium billiard the number of bouncing ball modes grows as N^^iE) ~ cE'^l'^ [ ]7D| , 



Therefore the bouncing ball mode contribution to S\{E,A) is, according to equation (^31) , 
proportional to and thus of the same order as the expected rate of quantum ergodicity 

for the quantum-ergodic eigenfunctions. One therefore expects for all domains in position space 
a rate of E~^l^. We have investigated the rate of quantum ergodicity for the stadium billiard 
using the small domains shown in fig. ^) and for larger domains shown in fig. |^ The results 
of the fits of S'^[E) = aE~^^^^^ to the numerical data for Si{E, XdJ are given in table 
Let us first discuss the rate for the small domains shown in fig. |^). For the domains Di and 
D2 which lie inside the rectangular part of the billiard the rate is in very good agreement with 
£;-i/4 g^^^ both for the domain D3, which lies on the border between the rectangular part and 



the quarter circle, and in particular for domain D4 which lies inside the quarter circle, one finds 
a slower rate than expected. This is a behavior which one would expect for a billiard with a 
much faster increasing number of bouncing ball modes. 



Figure 6: Plot of Si{E,xd,) for different domains Di for the stadium billiard using 
the first 6000 eigenfunctions, see fig. ^) for the location of the domains Di. The 
inset shows the same curves in double-logarithmic representation together with a 



fit of eq. (H) 



We see three possible explanations for this behavior of the rate for the stadium billiard. First, 
the counting function A'^bb(-E) for the bouncing ball modes might increase with a larger exponent 
than 3/4, A^bb(-E') ~ cE'^, [3 > 3/4. This would contradict the results in |7^, derived by 
independent methods. Moreover, the exponent j3 was tested numerically in up to energy 
E ^ 10000 and we found very good agreement with /3 = 3/4. Even if we relaxed the criteria for 
the selection of the bouncing ball modes drastically, the exponent did not change significantly, 
only the prefactor c increased. Therefore we think that this first possibility is clearly ruled out. 
Secondly, the rate for the quantum-ergodic eigenfunctions might not be proportional to E~^^^, 
but has a slower decay rate. Then we have to assume a position dependence of the rate, in 
order to explain the different behavior for the different domains: in the rectangular part of 
the billiard the rate has to be proportional to E~^^^ to explain the value of e obtained for 
the domains Di and D2. Whereas inside the quarter circle the rate of decay has to decrease 
as S[{E,xd4,) ~ i''{A)E^^'^^, in order to explain the value of e obtained for and D4. A 
priori such a dependence of the rate of the quantum-ergodic eigenfunctions on the location of 
the domain in the billiard is not impossible. If this is the case then one should observe no 
dependence of the rate on the volume of the domain D, as long as one stays in the same region 
of the billiard. E.g. the rate for a domain like Dq, which contains Di and D2 and is far enough 



domain 


rel. area 


e 


a 




biA) 


1 


0.015 


+0.009 


0.041 


0.0539 


0.0000 


2 


0.015 


+0.012 


0.041 


0.0564 


0.0000 


3 


0.015 


+0.033 


0.035 


0.0533 


0.0008 


4 


0.015 


+0.095 


0.029 


0.0492 


0.0047 


5 


0.015 


+0.020 


0.039 


0.0551 


0.0004 


6 


0.278 


+0.070 


0.137 


0.1401 


0.0233 


7 


0.433 


+0.111 


0.118 


0.1071 


0.0395 


8 


0.557 


+0.168 


0.089 


0.0292 


0.0634 


9 


0.696 


+0.188 


0.098 


0.0384 


0.0827 


10 


0.681 


+0.084 


0.176 


0.2474 


0.0295 



Table 2: Rate of quantum ergodicity for the stadium billiard with domains Di as 
shown in figs. 0and|^. Shown are the results for e and a of the fit Sf^{E) = aE~^^^~^^ 
to Si{E, XdJ- Also tabulated are the values for the relative area of the corresponding 
domains, and the results i^'ixoj and b{A) of the fit of the model ( ^31) to Si{E, XdJ- 



Figure 7: Domains in the a = 1.8 stadium billiard used to decide between the 
different explanations for the slow rates in the stadium billiard. 



away from the quarter circle, should be the same as the one for Di and D2. 
The third possible explanation for the observed behavior of the rate is that there exist more 
not quantum-ergodic eigenfunctions which have a larger probability density in the rectangular 
part than in the quarter-circle, and which are not bouncing ball modes. Alternatively the 
reason could be a subsequence of density zero of quantum-ergodic eigenfunctions, which has a 
sufficiently increasing counting function and a slow rate, see the remark at the end of section 
[^.4| . In both cases the model for Si{E,A) discussed in section O, which we already used in 
the case of the cosine billiard, would be applicable. In contrast to the second possibility in this 
scenario one expects a dependence of the rate of Si{E, xd) on the volume of the domain D, as 
in the case of the cosine billiard. 

To decide which explanation is the correct one we studied the rate for a number of large 
domains shown in fig. |^. With these domains one necessarily comes closer to the boundary of 
the billiard. To rule out the possibility that the observed behavior of the rate is due to the 
infiuence of the boundary, and not due to the dependence on the volume and location of the 
domains, we computed in addition Si{E,xd) for the small domain which is close to the 
boundary. 

The results are also given in table |^ and some examples of Si{E, XOi) for these large domains 
are shown in fig. |^. As for the cosine billiard, we also found that for large domains at small 
energies the rate may be much faster than at higher energies which is nicely seen in fig. § for 
the domains and D^. This effect is due to the influence of the boundary, as we will discuss 
in section |3.1.4| ; here we only note that the boundary influence vanishes for large energies. 
The observed rate of quantum ergodicity displays a strong dependence on the volume of the 
domain D, whereas the location, as long as one stays inside the rectangular part, has no 
influence. E.g. for the domain which contains Di and D2, one gets a much slower rate than 
for Di and D2- In contrast to the rate for the small domain near the boundary is rather 
close to the one for Di and D2- The slightly slower rate for is due to the smaller energy 
range for which we have computed Si{E^xd^)- A flt of Sl^{E) = aE'^^"^^^ to Si{E,xdi) and 
Si{E, XD2) using the flrst 2000 eigenfunctions gives an e of 0.022 for Di and 0.011 for D2, which 
is of the same magnitude as the result for D^. Moreover the rate decreases monotonically with 
increasing area of the domains D^, as long as they are inside the rectangular part R of the 
billiard. 

The domain Diq is interesting because it extends over both parts of the billiard. The enhanced 
probability density of the exceptional eigenfunctions in the rectangular part is partially compen- 
sated by the lower probability density in the quarter circle. Therefore one expects a rate similar 
to a domain in the rectangular part with relative area (vol(Dio) — 2vol(i^ion(n\i?)))/ vol(f2) = 
0.371 .... This relative area lies between the values for Dq and Dj, and indeed the rate for Diq 
lies between the rate for Dq and 1)7 too. 

These results strongly support the third explanation, i.e. the existence of a large density zero 
subsequence which is responsible for the deviations of the rate from E~^^'^. "Large" means 
that the counting function increases sufficiently strong to cause the rate to deviate from the 
expected behavior. 

To test this conjecture quantitatively one has to compare the numerical data with the conjee- 



Figure 8: Plot of Si{E,xd) for large domains (see fig. ^ for the a = 1.8 stadium 
billiard using the first 2000 eigenfunctions. The inset shows the same curves in 
double-logarithmic representation together with a fit of eq. (^Sf). For the domains 
Dj and in particular for domain a sharp transition from a fast to a slower decay 
of the rate is visible. This effect is due to the boundary and will be explained in 



sec. 3.1.4. 



tured behavior 



^modei^^^ A)= [I- cE'P) u\A)E-^l^ + h{A)E-^ . (52) 



Since this model contains the four free parameters c, /5, v'iA) and h{A)^ the numerical fit is 
not very stable. Therefore it is desirable to get some additional information from a different 
source. 

To this end we plotted di{n) = (■?/'„, XOii^n) —^oli^D^)/ vol(r2) for domain 1^4 which shows a slow 
rate, see fig. ^ and divided the spectrum into two parts by inserting a curve —caE~^^^^ and a 
curve c^E~^/^. The part of the spectrum between the two curves corresponds to the quantum- 
ergodic eigenfunctions with the optimal rate ~ and the part above and below the curves 

corresponds to the not quantum-ergodic eigenfunctions or to quantum-ergodic eigenfunctions 
with a slower rate than E~^^^. By computing the counting functions for these two subsequences 
we get a further criterion for distinguishing between the two possible scenarios for the behavior 
of the eigenfunctions discussed above. If the rate of quantum ergodicity for all quantum-ergodic 
eigenfunctions is slower than inside the quarter circle, then the fraction of eigenfunctions 

which lie below or above the two curves should grow proportional to E. If the deviation of 
the rate is due to a not quantum-ergodic subsequence of density zero, or a quantum-ergodic 
subsequence of density zero with exceptionally slow rate, then the number of states which lie 
below or above the two curves should grow like N"{E), i.e. slower than E. 
Proceeding in the described way, we find that the majority of the exceptional states have values 
of d4^{n) which lies below the lower curve —ca_E~^/^. A numerical fit for their counting function 
gives N"{E) = 0.06 E^'^^. The exponent is very stable under slight variations of the constant Cd 
which determines the curve. Up to E ^ 10000 corresponding to the 2000th state, the counting 
function even has an almost linear behavior. The nature of these states will be discussed below. 
The numerical result N"{E) = 0.06 E'^-^^ allows to determine the parameters c = 0.06 ^^'|^^^ = 
0.29 ... and /5 = -0.07 in the model (H) giving 



^modci^^^ A) = (l- 0.29 u'{A)E-'/^ + b{A)E-''-''^ . (53) 

We have now eliminated two of the four free parameters, and can therefore test this formula 
effectively with the numerical data. The results for v'{xDi) and h{xDi) are also shown in table 
0, and for three large domains the plot of Si{E, xd) and the corresponding fit S'™°'^'^'(£', xd) is 
shown in fig. |^. 

The agreement of the fits with the numerical data is very good. Moreover the values for v'{Di) 
and h{Di) are very reasonable: For ^'{Di) one expects that it depends on the volume of Di only, 
and not on the location. This is very well confirmed for the domains Di-D^, which have the 



same volume, where i''{Di) stays almost constant. According to the results in [0 one expects 
that v'{Di) increases with increasing volume of Dj, for small vol(Di), then reaches a plateau, 
and then finally decreases for very large domains. This behavior is not observed, the values 



for v'{Di) rather oscillate. The most striking difference occurs between u^Dg) and u^Dio), 
because they have approximately the same volume. We furthermore find that the behavior of 
v'{Di) is completely analogous to that of a,. The behavior of h{Di) is in perfect accordance with 
what one expects for a sum of quantum limits which are concentrated on the rectangular part 
of the billiard. The values increase when moving Di into the quarter circle, and they increase 
with increasing volume of Di, as long as Di lies entirely inside the rectangular part. For Dio 
the parameter 6(-Dio) takes an intermediate value between hi^Do) and b{Df), as one expects. 
The inclusion of the factor (1 — 0.29 in eq. (|53| ) turned out to be necessary to get 

satisfactory results. The contribution of E~^'^'^ cannot be neglected in the present energy range 
because of the small exponent. Without this factor we obtained for some of the domains 
negative values for v'{Di), which is impossible because S[{E,A) is by definition positive. 
This also sheds some light on the limitations of such a simple model like (0). For the exponent 
(5 = 0.07 only the order of magnitude is known for sure, the constant c from N"{E) might still 
vary, and nothing is known about the behavior of the higher order contributions to S'l and S". 
In view of this it is surprising how good this model fits with the numerical data. We believe 
that this gives a very strong support for the underlying conjectures, namely that a density 
one subsequence of quantum-ergodic eigenf unctions has a rate S[{E,A) ~ cE~^^^, and the 
deviations in the rate of Si{E, A) from this behavior are due to a subsequence of density zero. 
As mentioned in section a behavior S2{E, A) ~ cE~^^'^ ^^^ voUn) ^i/2j ^j^^ stadium billiard 
is claimed in |T^. We have tested this both for the small domains Di and D2, which are not 
influenced by the bouncing ball modes and also for some larger domains. However, the resulting 
fits clearly show that this result does not apply to our numerical data, see fig. H We also tested if 



this result applies to the quantum-ergodic subsequence, i.e. S[{E, A) ~ cE y ( -K^/^) , 



by replacing the term E in equation (|53|) by E ^/'^^hi{^ voi^n) ^ 1/2 j _ Again we find that 
from our numerical data that this possibility is excluded, at least for the energy range under 
consideration. For the stadium billiard it is known that the asymptotic behavior of the classical 
autocorrelation C(r) ~ 1/r, which leads to S2{E,A) ~ cE~^/'^ ln( ™^^^'' -E-'^/^) according to |jl6 



sets in rather late. So it would be very interesting to compare the results with those obtained 
by inserting the numerically computed autocorrelation function in the integral in p^). 
We now return to the question of what type these additional subsequences of eigenfunctions 
are. As additional information for the model, the counting function for the number of states 
for which (t/^^, XDii^n) — ™oiffif smaller than —c^E~^/'^ has been used. For comparison we 
have carried out the same procedure for the observable 1 — which corresponds to the 
complete quarter circle. As expected the bouncing ball modes appeared in both subsequences, 
but additionally a considerable number of other types of eigenfunctions showed up. In fig. [TO we 



Figure 9: Plot of S2{E,XDi) for the domains Di and D2 in the stadium billiard. 
The dashed lines show the fit of the conjectured behavior cE'^/"^ In^ ™'^^^-* ^'^/^) to 
S2{E,xdJ- The result of the fit shows that the numerical data for the first 6000 
expectation values cannot be described with this rate. 

Figure 10: Four examples of the exceptional eigenfunctions showing localization in 
the rectangular part of the stadium billiard, which are not bouncing ball modes, a) 
n = 1643, h) n= 1652, c) n = 1797 and d) n = 1834. 



show some examples of such eigenfunctions. They all show a reduced probability density inside 
the quarter circle, but their structure is essentially different from the bouncing ball modes. 
Their semiclassical origin are maybe periodic orbits bouncing up and down between the two 
perpendicular walls for a long time but then leaving the neighborhood of the bouncing ball 
orbits in phase space. At least it seems difficult to associate short unstable periodic orbits to 
the patterns in the shown states, because the lines of enhanced probability do not always obey 
the laws of reflection, or they look too irregular. 

A further test of the hypothesis that a density zero subsequence is responsible for the slow rate 
is provided by varying the length a of the billiard. Here we used the first 2000 eigenfunctions 
for both the a = 0.5 and the a = 4.0 stadium billiard in addition to the results for the a = 1.8 
stadium based on 6000 eigenfunctions. We have chosen three different domains for these three 
systems: domain A lies within the rectangular part of the billiard, domain B is centered at 
X = a and domain C is located in the quarter circle. The results for the rate of quantum 
ergodicity are shown in table For different parameters the quantities b{Di) change, and 
therefore the weights of the different contributions to Si{E,A) in equation (^31). For smaller 
a the relative fraction of the volume of the rectangular part, vol(-R)/ vol(f2) becomes smaller. 
Therefore one expects that for smaller a the influence of the not quantum-ergodic subsequences 
to Si{E,xd) becomes stronger in the rectangular part, and weaker in the quarter circle. This 
is nicely seen in the numerically found behavior of the rate for the domains A and C shown in 
table ^, which confirms our hypothesis. 

To summarize our results for the stadium billiard, we have shown the existence of a large, but 
density zero, subsequence of eigenfunctions which have an enhanced probability distribution 
on the rectangular part of the billiard but having a different structure than the bouncing ball 
modes. We demonstrated that the observed effects are due to the influence of this subsequence 
of density zero. This subsequence shows a different behavior than the majority of quantum- 
ergodic eigenfunctions for which our results imply a uniform rate of E~^^^. Clearly we cannot 
decide if this exceptional subsequence will ultimately be not quantum-ergodic, or if it is a 
quantum-ergodic subsequence with a exceptional behavior of the rate. We can only say that 
on the presently studied energy range up to £J ~ 30000, i.e. up to the 6000th eigenfunction, 
they behave not quantum-ergodic. 



system 


domain A 


domain B 


domain C 


stadium (a = 0.5) 


+0.111 


+0.062 


+0.056 


stadium (a = 1.8) 


+0.009 


+0.033 


+0.095 


stadium (a = 4.0) 


-0.008 


+0.031 


+0.095 



Table 3: Results for e of the fit of Sf^{E) = aE^^^^'^^ to the numerically obtained 
Si{E,XDi)i for stadium billiards with different parameter a for three different do- 
mains A, B and C . Domain A lies within the rectangular part of the billiard, domain 
B is centered at x = a and domain C is located in the quarter circle. 



3.1.3 Cardioid billiard 



The cardioid billiard is probably the most "generic" one of our three billiards, in the sense that 
it possesses no two dimensional family of periodic orbits like the bouncing ball orbits. One 
might therefore expect a priori a better rate of quantum ergodicity than for the other billiards. 
We have computed Si{E,XDi) for five small domains, see fig. ^c), by using the first 6000 
eigenfunctions up to energy E ^ 32000, and for three larger domains, see fig. by using the 
first 2000 eigenfunctions. The results are displayed in figs. 11 and 12 . To determine the rate a 
fit of Sf^{E) = aE~^f^^^ has been performed, and the resulting values for a and e are listed in 
table |. 

We find that domain gives the lowest rate of quantum ergodicity for the small domains Di- 
. This is caused by a considerable number of eigenfunctions showing an enhanced probability 
as in fig. ^) along the vertical orbit AB. For domains Di, D2 we also find a slower rate than 
for the other regions 1^4, D5; in this case the slower rate seemingly cannot be attributed to one 
type of localized eigenfunctions. 



Figure 11: Plot of Si{E, XOi) for different domains Di for the cardioid billiard using 
the first 6000 eigenfunctions, see fig. ^) for the location of the domains Di. The 
inset shows the same curves in double-logarithmic representation together with a 
fit of Sf'{E) = aE"l/4+^ eq. 

The larger domains show a slower rate than the small domains, but the rate is not monotonically 
decreasing with the area of the domain. The rate for the largest domain Dg is even of the same 
order of magnitude than the one for D3, especially if one takes the smaller energy range for Dg 
into account. This slower rate is probably caused by the existence of different not quantum- 
ergodic subsequences with quantum limits fik in different regions of the billiard. For each of 
the domains the influence of these subsequences is different and therefore one observes different 
rates. 

A quantitative test in a similar way as for the other billiards using a model for Si{E,A) is 



domain 


rel. area 


e 


a 


1 


0.01722 


+0.047 


0.028 


2 


0.01722 


+0.039 


0.037 


3 


0.01722 


+0.064 


0.046 


4 


0.01722 


+0.007 


0.048 


5 


0.01722 


+0.009 


0.042 


6 


0.18674 


+0.098 


0.125 


7 


0.33104 


+0.115 


0.140 


8 


0.50930 


+0.071 


0.213 



Table 4: Rate of quantum ergodicity obtained from a fit of Sf^{E) = aE ^/'^+^ to 
Si{E, XOi) for the cardioid billiard with domains Di as shown in figs. |^) and |12. 



Figure 12: Plot of Si{E,XDi) for larger domains for the cardioid billiard using the 
first 2000 eigenfunctions. Also shown are fits to eq. (|48|) for the corresponding energy 
regions. 



very difficult, because the deviations from the conjectured optimal rate is not only due to one 
subsequence. But the results for and clearly shows that here as well one has a density 
one subsequence of quantum-ergodic eigenfunctions with rate S[{E,xd) ~ ^' 
hope to return to the problem of determining the not quantum-ergodic subsequences and their 
quantum limits in the future. 

The cardioid billiard is the only system we have studied to which the result (|33D should be 
applicable. But for most of the domains the rate is much slower than the predicted one. Only 
the domains 4 and 5 show the expected rate. Therefore we have computed for these domains the 
factor p{A) in eq. (^). For the computation of p{A) the variance of {xDi)i — ^"oiln) ^ function 



vol(n) 

of / has been computed using trajectory segments of length I of a generic trajectory {qit)}. The 
quantity {xDi)i = j Jo XOiiqit)) dt is the relative length of the trajectory segment lying in the 
domain A- By ergodicity we have lim/^oo(XD,)/ = ^py- The variance of {xd,)i - ^^fp^ 
decreases like p{A)l~^. 

Using the corresponding results in equation (p3D we obtain S2{E,xd4) = 0.0062 and 
S2{E,xd5) = 0.0074 .E"^/^. These numbers have to be compared with the result of a fit 
SfiE, A) to S2iE,XD J- We ohtainS f{E,XD,) = 0.0036 ^-O-^^ and xds) = 0.0031 ^-o-^^ 

On sees that the theoretical prediction is too large by a factor of approximately 2. This de- 
viation might be related to the factor g in (p3D, which counts the mean multiplicities in the 
classical length spectrum. In the cardioid billiard the asymptotic value g = 2 is reached very 
late, for the shorter periods one rather has g ^ 1, which would lead to a better agreement of 
eq. ( p3D with the data for D4 and D^. 

For a better understanding it seems necessary to check in detail, whether any of the assumptions 
leading to eq. (|33|) is not fulfilled for the domains of the cardioid billiard. It would also be very 
interesting to investigate if the slower rates can be described using the expression in terms of 
the classical correlation function. We will leave these questions for a separate study. 

3.1.4 The influence of the boundary 



Figure 13: In a) we show a three dimensional plot of the sum \E'£'(a;, = 
Iv^ J2e„<e \i'i^yy)\'^ involving the first 250 eigenfunctions of the a = 1.8 stadium 
with odd-odd symmetry. The pictures on the right show a cross section \E'£;(1,?/) 
for using the first b) 250 and c) 1000 eigenfunctions. The dashed curves in b) and 
c) display the evaluation using the first two terms in formula (|5^. These results 
are used to explain the fast rate in the low energy range for the stadium billiard for 
large domains. 

In all three billiards we observe the phenomenon, that for large domains Si{E, xd) decays faster 
at low energies than at high energies. This can be seen in fig. |^ for domain Dg in the cosine 



billiard, in fig. § for domains and Dg in the stadium billiard and in fig. for domains 
Dq, D-j and in the cardioid billiard. The other large domains we studied showed the same 
behavior. The only exceptions are the domains in the cosine billiard and in the stadium 



billiard, which consist of the whole rectangular part. For these domains no faster rate at low 
energies is visible. 

Qualitatively this behavior can be understood by the vanishing of the probability density 
of the eigenstates at the boundary due to the Dirichlet boundary conditions. Be- 
cause of the normalization of ipniQ) the reduced probability density at the boundary has to be 
compensated by an enhancement of the probability density inside the billiard, which leads to 
larger oscillations of the probability density near the boundary. 

Let us assume that this compensation of the probability density takes place in a strip along 
the boundary of a few de Broglie wavelength width. Then the integral of the probability 
density \ipn{(l)\'^ over a domain D feels the influence of the boundary only up to a certain 
energy, proportional to the inverse square of the distance between D and the boundary dQ. 
Furthermore the boundary influence will be proportional to the overlap of D and the strip at 
the boundary. This overlap decreases like and therefore Si{E,xd) should decrease 

with such a rate at low energies. So the assumption that the compensation takes place in a 
small strip along the boundary leads exactly to the behavior we observe. Moreover a domain 
like Dg which extends to the boundary dQ will not feel any influence, because the boundary 
effect is compensated entirely inside this domain. 

To justify our assumption on the range of the boundary influence we refer to the following 
result on the asymptotic behavior of the summed probability densities on a two-dimensional 



Riemannian manifold with C°°-boundary |^0|, Theorem 17.5.10] 



J2\Mlf-^E-- > VE + R(,,E), (54) 



En<E 



where d{q) is the shortest distance of the point g G to the boundary. The remainder -R(g, E) 
satisfies the estimate |-R(g, E)\ < C\fE. The second term in describes the influence of the 



boundary, for d{q) — > the term tends to —E /{An) and cancels the contribution from the first 
term, such that the boundary conditions are fulfilled. In fig. 0a) the normalized sum 

'^E{x,y) = \Mx,y)\^ , (55) 

^ ^ E„<E 

is displayed for the stadium billiard, using the first 250 eigenf unctions. One nicely sees how 
the probability density is forced to vanish at the boundary, and how the compensation leads to 
large oscillations near the boundary. In fig. |T3| b) and c) we show two cross sections through 
the function (|55|) at two different energies, and compared it to the result one gets from the first 
two terms on the right hand side of (0). The agreement is quite impressive, especially near the 
boundary (y = and y = 1). So although the stadium billiard does not have C°°-boundary, 
the result ( |5^ seems to remain valid. One furthermore observes that with higher energies the 
y-range on which the agreement is excellent increases. 

The averaged probability density (^) shows exactly the behavior we assumed for the individual 
wavef unctions, in order to explain the fast rate of quantum ergodicity at low energies for 
domains near the boundary. The influence of the Dirichlet boundary condition is concentrated 
near the boundary, and it decays at a length scale proportional to the de Broglie wavelength. 
So with the help of eq. (0) one gets a good qualitative understanding of the boundary influence 
on the rate of quantum ergodicity. 



In order to try to get a quantitative understanding we used eq. (Q) to derive as in a mean 
eigenfunction which incorporates the boundary influence 



\^n{q)? - ^^^^^^^^^ (l - Jo(2rf(g)v^)) . (56) 



2\/ En 



Integrating this expression over a domain D should give for the expectation values the mean 
value plus the corrections due to the boundary of xd- By incorporating this into Si{E, xd) one 
obtains an expression, which we compared with our numerical data. Although ( ^61) implies a 
faster decay rate at low energies, it is not as strong as the numerically observed one. This devi- 
ation must be caused by considerable fluctuations of the boundary influence on the individual 
states ipn around the mean influence described by (0) and (|56|) . 



3.2 Quantum ergodicity in momentum space 

Up to here we have investigated the behavior of the wavefunctions in position space only. Now 
we turn our attention to the rate of quantum ergodicity in momentum space, which is studied 
here for the first time numerically. 

Quantum ergodicity predicts that the angular distribution of the momentum probability dis- 
tribution |'?/'n(p)P tends to l/(27r) in the weak sense, see eq. (p8|). Therefore we study an 
observable with symbol Xc(e, Ae) ip) whose expectation value gives the probability of finding the 
particle with momentum-direction in the interval ]6 — A6'/2, 6 + A6'/2[. Recall that Xc(0,Ae)(p) 
denotes the characteristic function of the circular sector C{9,A9) = {p G | a.Tcta.n{py/px) G 
]9 — A9/2, 9 + A9/2[}, and the classical mean value of Xc(e,Ae)(p) is A9/ {2n). 
Only eigenfunctions of odd parity of the not desymmetrized systems are considered here due to 
our method of computing the Fourier transformation directly from the normal derivative m„(co') 
of the eigenfunction ipnil)- From Green's theorem one easily finds the formula 



Mp) = -TT^^ ^ / e-'''(-)P«„(u;) do. , (57) 



where q{oj) denotes a point on the boundary dVL. The advantage of this formula is that it 
allows to compute the Fourier transform directly from Un{oj), which can be obtained using the 
boundary integral method. For desymmetrized systems, like the ones considered here, one uses 
an appropriate Greens function which vanishes at the lines of symmetry, and therefore removes 



them from the boundary integral, see e.g. |j71[. This reduces the computational effort, but one 
does not get the normal derivatives on these parts of the boundary of the desymmetrized system. 
Therefore our results for the rate of quantum ergodicity in momentum space are sufficient to 
rule out the possibility of a totally different behavior in momentum space than in position 
space. Since the rate for all eigenfunctions cannot be faster than the one for a subsequence of 
positive density, we get a lower bound for the rate of the full system. 

The time reversal invariance leads for the Fourier transformed eigenfunctions to the symmetry 

ipn{—p) = '4'nip)- Therefore \ipn{—p)\^ = |V'n(p)P, and this reduces the angle interval we have to 
study to [0,7r[. The additional reflection symmetries in the considered billiards further reduce 
the relevant angle interval to [0,7r/2[. 

For our numerical computations we have chosen five equidistant intervals, centered at 9i = 
{i - 1/2) with i = 1, . . . , 5 of width = f^. 



system 


domain 


e 




1 


0.15 




2 


0.12 


stadium 


3 


0.15 




4 


0.09 




5 


0.18 




1 


0.050 




2 


0.075 


cardioid 


3 


0.026 




4 


0.079 




5 


0.076 



Table 5: Rate of quantum ergodicity obtained from a fit of Si^{E) = aE^^/^^^ to 
tlie numerically obtained function Si{E, Op[xciei,A0)]) ^^e different systems and 
angle sectors C{9i,A9). 



As in the case of quantum ergodicity in coordinate space, see eq. ( ^61) and fig. y, one observes 
large fluctuations of {4'n,Xc{9,Ad)'^n) — A^^/(27r) around 0. Therefore we again consider the 
cumulative version (|3T|) of the quantum ergodicity theorem, which reads in this case 



Si{E,Op[xc{eAe)]) 



NiE) 



A0 



\A{p)\^ dp - — 



C(0,A6») 



for E ^ oo 



(58) 



Figure 14: Plot of ^i(E, Op[xc(e„Ae)]) for 6^ = {i - 1/2) with i = 1, ... ,5 and 
A6 = for the stadium billiard using the first 2000 eigenf unctions. 

Figure 15: Plot of Si{E,Op[xc(e,,A0)]) ^^ = - 1/2)y^ with i = 1, ... ,5 and 
= for the cardioid billiard using the first 2000 eigenfunctions. 



The results for Si{E, Op[xc(6»,,Ae)]) are shown in fig. |T4| for the stadium billiard and in fig. [15 
for the cardioid billiard. In each case 2000 eigenfunctions have been used. For the cardioid 
billiard the inset shows a double logarithmic representation together with the fits of Sf^{E), 
eq. (^). For the cosine billiard no computations of the rate Si{E, Op[xc(6'i,Ae)]) momentum 
space have been performed. 

As in position space one expects that the rate is strongly influenced by not quantum-ergodic 
subsequences of eigenfunctions. For the bouncing ball modes in the stadium billiard one has 

lim / |Vn"(p)|Mj9= r Ae,9 + Ae[ 

K"^^ J ^ ll for ^e]9-Ae,e + A9[' ^ ' 




and so the coefficient t'bb(Xc{ei,Ae)) ^^e model (|43|) for Si{E,xc{9iAd)) given by 

The results for the rate of quantum ergodicity, characterized by e, are listed in table It 
turns out that the rate is slower than the rate of quantum ergodicity for the small domains 
in configuration space. Moreover the agreement of Si{E,Op[xc{9i,A9)\) with the fit is not as 
good as in the case of Si{E, Xd), in particular the fluctuations of Si{E, Op[xc{ei,A9)]) much 
larger than in position space. 

In the stadium billiard the interval 5, which corresponds the direction of the bouncing ball 
orbits, shows the slowest rate. But as we already noted in the discussion of the rate in position 
space, the bouncing ball modes alone cannot cause such a slow rate, because their counting 
function increases only as E^^^. So a considerable number of the additional not quantum- 
ergodic states which are responsible for the slow rate in position space must also have an 
enhanced momentum density around 7r/2. But the slow rates for the other angular intervals 
indicate that not all not quantum-ergodic states show this behavior in momentum space. 
For both billiards one observes that the order of magnitude of e in momentum space is the 
same as for the large domains in position space. Therefore the results are compatible with the 
results in position space, but the large fluctuations indicate that one has to go higher in the 
energy in momentum space than in position space. 



3.3 Fluctuations of expectation values 

Another aspect of great interest is how the expectation values {ipn, fluctuate around their 
mean value o"(A). Since the mean fluctuations decrease for large n, one has to consider the 
distribution of 



5. = ^•''•M>„)-<r(A) (g^j 
S2(E„,A) 

Here 5*2 (-E, A) = S 5*2 (-E", A) with S being a correction necessary to ensure that the distribution 
of has unit variance; see below for an explanation. So the question is whether a limit 
distribution P{^) of ^„ exists in the weak sense, i.e. 



N °? 

i^L^E^(^")= / aiOPiOd^ , (62) 

n=l 



where g{^) is a bounded continuous function. It is natural to conjecture that this distribution 
tends to a Gaussian normal distribution, 

P(0 = -^exp(-eV2) , (63) 
V /vr 



as in random matrix theory section VII]. Note that this is a conjecture for every observ- 
able, i.e. the asymptotic distribution should be independent of the special observable under 
investigation. For hyperbolic surfaces a study of -P(0 for an observable in position space is 



Figure 16: Cumulative distribution of ^„ = {{"ipn, Aipn) — cr{A.)) / y S2{En, A) for 
the stadium billiard for domain 4, A = XD4, and for the cardioid billiard with 
observable A = XD4, — XDs- In both cases we haven chosen n G [2000,6000]. The 
dashed curve corresponds to the cumulative normal distribution. The insets show 
the distribution of ^„ together with the normal distribution with zero mean and unit 
variance, eq. (|63|) (dashed curve). 



Figure 17: Cumulative distribution of ^„ = ((V'n, Atpn) — y S2{En, A) for the 

cardioid billiard, for the observable Xc{e,AS){p) in momentum space with 9 = 57r/20 
and A6 = vr/lO. The dashed curve corresponds to the cumulative normal distribu- 
tion. The insets show the distribution of C,n together with the normal distribution 
with zero mean and unit variance, eq. (|63|) (dashed curve). 



contained in [|T^], where a good agreement with a Gaussian normal distribution was observed. 



In |]T6| -P(0 was studied for the Baker's map and the hydrogen atom in a strong magnetic field, 
and a fair agreement with a Gaussian was found. 

However, already from the plots of di{n) shown in fig. |]it is clear that the fluctuations are not 
symmetrically distributed around zero, but have more peaks with large positive values. The 
reason is that di{n) = {ipn, XOii'n) — ^°oi(n) to satisfy the inequality 



vol(A) ^ I . vol(A) ^ ^ vol(A) 



vol(^]) - vol(^]) - vol(^]) 

This already indicates that the approach to an asymptotic Gaussian behavior could be rather 
slow. Therefore we have tested additionally for the cardioid billiard the observable A = XD4 — 
XD5 where the expectation values fluctuate symmetrically around zero, and one expects a faster 
approach to a Gaussian behavior. In fig. |16|a) we show the cumulative distribution 

/^(0 = ^#{n I a<0 (65) 

for domain D4 of the stadium billiard and in fig. |l6|b) In{0 is shown for the observable A = 
XDi — Xds in case of the cardioid billiard. In both cases all values of ^„ with n G [2000, 6000] 
have been taken into account, giving = 4000. For the rate S2{E,xd) we used the result of 
a fit to 5*1* (-E) = aE"'. The insets show the corresponding distributions of C,n in comparison 
with the normal distribution, eq. (|63|). Notice, that no further fit of the mean or the variance 
of the Gaussian has been made. Example a) is the case for which we have found the worst 
agreement with a Gaussian (of all the small domains we have tested). The observable chosen 
for b) gives very good agreement with the Gaussian distribution. In the case of XD4 in the 
stadium billiard there is a significant peak around ^ = —2, which is due to the bouncing ball 
modes, for which {4'n",XD44'n") is approximately zero, see fig. ^. Therefore one has a larger 
fraction with negative ^n"- For the distribution in case of the observable A = XD4 ~ Xdb of the 
cardioid billiard we obtain a significance niveau of 23% for the Kolmogorov-Smirnov test (see 
e.g. |7^) with respect to the cumulative normal distribution. 

We also studied the distribution of ^„ for the observables a{p,q) = a{p) = Xc{e,AS)i.P) i^ 
momentum space. For the stadium billiard the computed distributions show in the considered 



energy range clear deviations from a Gaussian, as one already expects from fig. |1^. The best 
result was obtained for the cardioid billiard for the interval given by i = 3 (with 9i = (i — 1/2)^ 
and A6 = j^) and is shown in fig. The agreement is quite good, the Kolmogorov-Smirnow 



test gives a significance niveau of 29%. 

There is one subtle point concerning the variance of the distribution of ^n- Since S2{E,A) 
does not represent a local variance around E, but a global one, it is necessary to take this into 
account to obtain for the fluctuations a variance of unity. If the rate behaves as 5*2 (i?. A) = aE"' 
then the correction S is given by S = a + 1, e.g. for a = —1/2 we have S = 1/2. See for 
a more detailed discussion on this point in the case of the distribution of the normalized mode 
fluctuations. 

Let us now discuss the influence of not quantum-ergodic sequences to the possible limit distri- 
bution. Assume that the rate for the quantum-ergodic states is S2{E,A) ~ aE~°' with some 
power a. If we have a subsequence of not quantum-ergodic states such that the total rate is 
S2{E,A) ~ a'E""" , we can have two different situations, either a = a', or a < a'. In the 
first case the not quantum-ergodic states have no influence on the limit distribution. In the 
second case where the not quantum-ergodic states dominate the rate, the normalization by a 
rate which is slower than the one of the quantum-ergodic subsequence implies that we have 

If one instead normalizes the fluctuations with the rate of the quantum-ergodic subsequence, 
S'2iE,A), 



a := , (66) 

S^(E„,A) 

with S2{En,A) = ESI^iE, A), then the limit distribution is determined only by the quantum- 
ergodic subsequence, independent of the behavior of the not quantum-ergodic subsequence. To 
see this we split (p2|) into the different parts 



1 N'{E) 1 r \ N"{E) 1 r ^ 

-WiE) ^ ^^^"^ - immE) ^ ^^^"'^ ^ immE) ^ ^^^-"^ • ^^^^ 

En<E E^i<.E E^//<E 

Since lim^^oo = 1, hniij^oo = and jj7^T.e„„<e 9{L") < ^^^^9(0^ one gets 

^ ^ E„<E ^ ^ E„,<E 

We conjecture that the fluctuations of the quantum-ergodic subsequence is Gaussian, and there- 
fore all fluctuations, when normalized with the rate of the quantum-ergodic subsequence, are 
Gaussian. 



4 Summary 

The aim of the present paper is to give a detailed study of the rate of quantum ergodicity 
in Euclidean billiards. We first have given a short introduction to the quantum ergodicity 
theorems in terms of pseudodifferential operators. We have shown that the quantum ergodicity 
theorems of Shnirelman, Zelditch, Colin de Verdiere and others are equivalent to a weak from 



of the semiclassical eigenfunction hypothesis for ergodic systems put forward in g |T0|, |Tl|, |T2l . 
That is, the quantum ergodicity theorem is equivalent to the statement that for ergodic systems 
the Wigner functions Wn{p, q) fulfill 

W„^. (p, q) ~ . 5{H{p, q) - , (69) 

for Erij oo and {rij} C N a subsequence of density one. 

Of great importance for the practical applicability of the quantum ergodicity theorem is the 
question, at which rate quantum mechanical expectation values {ipn, Aipn) tend to their mean 
value cr{A). Different arguments were presented previously in favor of an expected rate of 
quantum ergodicity Si(E,A) = 0(£'~^/'^+^), for all e > 0, in the case of strongly chaotic 
systems. In section we discussed the influence of not quantum-ergodic subsequences to 
the rate. If their counting function increases sufficiently fast, they can dominate the behavior 
of Si{E,A) asymptotically. Together with results from for the number of bouncing ball 
modes in certain billiards it follows that one can find for arbitrary 5 > an ergodic billiard 
for which Si{E,A) = 0{E^^). That is, the quantum ergodicity theorem gives a sharp bound, 
which cannot be improved without additional assumptions on the system. 
We furthermore developed a simple but powerful model for the behavior of Si{E,A) in the 
presence of not quantum-ergodic eigenfunctions, whose main ingredient is that the quantum- 
ergodic eigenfunctions should obey the optimal rate E~^^'^. The discussion shows that the total 
rate of quantum ergodicity can be strongly influenced by those exceptional subsequences. Not 
only that they can cause the rate to be much slower than E~^^*, they can lead as well to a 
grossly different behavior of 5*1 {E, A) at low and intermediate and at high energies. 
The numerical investigations are carried out for three types of Euclidean billiards, the stadium 
billiard (with different parameters), the cosine and the cardioid billiard. The results are based 
on 2000 eigenfunctions for the cosine billiard, and up to 6000 eigenfunctions for the stadium and 
the cardioid billiard. As observables we have used characteristic functions of different domains 
in position space and also a class of observables in momentum space, considered here for the 
flrst time. 

It turns out that the rate of quantum ergodicity in position space is in good agreement with a 
power law decay, Si{E,A) ~ The difference e between the exponent and 1/4 is found 

to be small for several domains and systems. However we also flnd a number of signiflcant 
examples showing a slow rate (i.e. e > and not small). These are discussed in detail and can 
be attributed to subsequences of locahzed eigenfunctions. 

For the cosine billiard we flnd that the rate agrees well with the expected rate, in particular for 
the small domains. However asymptotically the rate has to obey Si{E,A) ~ because 
the counting function of the bouncing ball modes increases as E^^^^. The asymptotic regime for 
the rate lies far beyond any presently computable number of energy levels. By incorporating 
the knowledge on the counting function obtained from our previous work, we tested our model 
(|43|) for the rate successfully for all the considered domains. 

For the stadium billiard the situation is more complicated: here the counting function of the 
bouncing ball modes increases as E^^'^ and therefore, as discussed in sec. |2.4| , cannot influence 



the rate. However, we flnd for the stadium billiard that the rate is for several domains in 
position space much slower than expected. After discussing and testing several possibilities 
our explanation for this observation is that in the stadium billiard there exist a much larger 
subsequence of eigenfunctions which have an enhanced probability density in the rectangular 
part of the billiard than just the bouncing ball modes. They nevertheless have density zero. 



but their counting function increases stronger than E^^'^. Of course, we cannot decide whether 
this subsequence either has a quantum hmit different from the Liouville measure, or if it is a 
quantum-ergodic subsequence, with an exceptionally slow rate. 

For the cardioid billiard we also have domains for which the rate is proportional to E~^^'^. But 
we also find significant exceptions, in particular for domain D3 the rate is much slower, and 
this can be attributed to a number of eigenstates which show localization along the unstable 



periodic orbit AB. For the cardioid billiard we also tested the result from |T6[, eq. (p3D, for the 
domains D4 and D^, for which the rate is closest to the optimal rate. However the semiclassical 
result does not agree with our numerical results for the rate. It would be interesting to study 
this in more detail. 

From our numerical results we obtain the following general picture: In the studied systems there 
is a quantum-ergodic subsequence of density one whose rate is S[{E,A) = 0{E~^^'^~^^). If one 
observes a slower rate of Si{E,A) by using all eigenf unctions, this is caused by a subsequence 
of density zero, whose counting function increases stronger than E^^'^. These exceptional eigen- 
functions show localization effects and and probably they tend to some non ergodic quantum 
limit. However we cannot rule out the possibility that they are quantum-ergodic but with a 
much slower rate than the majority of eigenfunctions. 

We have found furthermore an effect due to the boundary conditions. For domains lying next 
to the boundary we observed that the rate may be considerably faster at low energies. The 
qualitative explanation of the phenomenon is that the probability density of the eigenstates 
show enhanced fluctuations near the boundary because of the boundary conditions. 
Using an observable depending only on the momentum, we studied quantum ergodicity in 
momentum space too, which is done here for the first time to our knowledge. We find that in 
general the rate of quantum ergodicity is of the same magnitude as for the large domains in 
position space. Furthermore the oscillations of Si{E,A) are larger in momentum space, which 
might indicate that one has to go higher in the energy in momentum space than in position 
space. 

We also studied the distribution of the suitably normalized fluctuations of {tpn, Aipn) —cr{A), see 
eq. (l6l|), both for operators in position space and in momentum space. For the observable A = 
XD4—XD5 in the case of the cardioid billiard we flnd very good agreement with a Gaussian normal 
distribution and in the case of observables depending only on the momentum good agreement 
is found. However for the stadium billiard (and also domain D3 for the cardioid billiard) we 
clearly flnd that again subsequences of not quantum-ergodic states may have a considerable 
influence. If they dominate S2{E,A), the distribution will tend to a delta function due to the 
normalization by a/ S2{E, A). But when normalizing instead with the rate of the quantum- 
ergodic states, S'2(-E, A), we expect a universal Gaussian distribution of the fluctuations. 
As possible investigations for the future it seems very interesting to study whether the expression 
given in |T^ for the rate in terms of the classical correlation function can describe our numerical 
results. In particular for the cardioid billiard a more detailed investigation along these lines 
seems promising as this system is the most "generic" one from the three studied systems and 
we flnd both the optimal rate and also clear deviations. The present paper also shows that a 
detailed understanding of the phenomenon of scarred eigenfunctions is necessary because these 
clearly affect the rate of quantum ergodicity. 
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Appendix 

A Kohn-Nirenberg quantization 

In the mathematical literature one often prefers a different quantization procedure, sometimes 
called the Kohn-Nirenberg quantization [^0|, ^ , and the literature on quantum ergodicity often 
refers to this quantization procedure. To the symbol a G 5"'"(]R^ x Q) one associates the operator 

Op™[a]/(g) := ^ J e'^«a(p, q)f{p) dp , (70) 

where f{p) := JqG~^'^^ f{q) dq is the Fourier transform of /. The principal symbol is defined in 
the same way as before, i.e., if a ~ ^^lo '^m-fc, then the leading term am is called the principal 
symbol, a™ (Op™ [a]) = a^. 

The usual quantum ergodicity theorem is now the same theorem as we have stated it, but with 
the Kohn-Nirenberg principal symbol instead of the principal symbol corresponding to the 
Weyl symbol which we have used. But it is well known, see pT] , pO| , that if a G ^"^(M" x M"), 
then the Weyl symbol of the Kohn-Nirenberg operator belongs to the same symbol space, 
WfOp^'^fa]] G S''^(M" X M"), and that the principal symbol coincides with the Kohn-Nirenberg 
principal symbol, 

a(0p™H)=a™(0p™[a]). (71) 
Therefore the two formulations of the quantum ergodicity theorem are equivalent. 



B Generalizations of the quantum ergodicity theorem 

Assume we have given a quantum limit fik on Si, that is we have a subsequence of eigenfunctions 
{i^n,}j&^, such that 

\im{ipn,,Aipn,) = / fikip,q)(riA){p,q)dfx , (72) 

J^OO J 

for all A G 5'°i(fi). We want to discuss the lift of from Si to the whole phase space. 
To this end we express the expectation values for an operator of arbitrary order m G M by 
the expectation values of an operator of order zero. This can be achieved by using the fact 
that for every m G M, (— A)^ is a pseudodifferential operator of order m with principal symbol 
a((-A)f) = {a{-A)p = if(p, g)?, see |2|, |2|. By multiplying an operator A G 5" (R™) 



of order m with the operator (—A) 2 ^ which is of order —m, we get an operator (—A) 2 A G 
^^(IR") of order zero. For the expectation values of A we therefore have 

m m 

AV-n,) = En] (V-n,, (-A)-^^^^^.) , (73) 

and on the right hand side we have an operator of order zero. 

The principal symbol of {—A)~^A is according to eq. (^TJ) given by a{{—A)^^)o'{A) = 
H{p, q)~~a{A), and since by definition H{p, g) = 1 on Ei we obtain from (|72D and (^) 

lim ' {^n, , A^n, ) = / fikip, q) cr{A) (p, q) d/x . (74) 

Si 

Thus eq. (|7^ provides the extension of the quantum ergodicity theorem to pseudo-differential 
operators of arbitrary order m. 



C Connection to the semiclassical eigenfunction hypoth- 
esis 



By introducing the definition of the Liouville measure /x, equation ([7^) can be written as 

{^n,,A^n,) jj /i,(p, q)a{A){p, q)^AE^^-A dpdg . (75) 

If one uses the homogeneity of cr(A), i.e. En^a{A){j), q) = a{A){EnjP, q), and performs a change 

_ 1. 

of the momentum coordinates from p to En/p one obtains 

_ 1 

vol(Si)E„^ 



where furthermore the homogeneity properties of H{p, q) and of the delta function have been 
used. In terms of the Wigner functions Wn^ corresponding to ipnj eq. (|7^) reads 

(y{A){p, q) Wn^{p, g) dpdg ~ [[ (r{A){p, q) fi^iE'^p, q) "^^^^^'^^ ^^^^ ^p^^ ^ (77) 

J J vol(Si)E„^^^ 

where (j{A){p^ q) can be any function homogeneous in p of degree m, for some arbitrary m G M. 
But since the set of all polynomials in p is already dense in C°° {E? x Vt) the set of homogeneous 
functions in p is dense in C°°(]R^ x fi), too. Therefore one gets 

u/ ( \ (Tp-'2 ^ S{H{p,q)-E ) 

Wn, [p, q) ~ MEn/p, q) . (78) 

vol(Si)E4. 



— — 1 

Note that vol(Ei)£'^^. = vo1(Ee„^. ), and if we extend q) from Ei to the whole phase space 

by requiring it to be homogeneous of degree zero, q) — l^kip/ ^/H{p, g), g), then we finally 
can write 

(P, q) ~ /^fcb, q) 77^ ^ for J ^ oo . (79) 

vo1(Ee„.) 

for a subsequence {rij} C N of density 1. This shows that the quantum ergodicity theorem is 
equivalent to the semiclassical eigenfunction hypothesis for ergodic systems for a subsequence 
of density one. 
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